%% Rutherford Scattering

clear all
close all

Z1=2;   %Charge on alpha
Z2=79;  %Charge on Gold

m_e=9.1*10^-31;     %kg
e=1.6*10^-19;       %C
k=8.9*10^9;         %kg m^3 s^-2 C^-2
A=k*Z1*Z2*e^2/m_e;  %m^3/s^2
V_0=2*10^7;         %m/s

% Setup the ODE solver
P=ODESolver();
P.f=@(X)CentralPotential(X,A)

% propogate past the nucleus
Z_0 = -1e-8;        % (m) start back where force is small
P.T=-2*Z_0/V_0;     % (s) go past nucleus
P.dt=P.T/100;       % (s) Divide into 100 time steps

%Max impact parameter
Max_b=1e-9;     %m

figure
plot(0,0,'yo','linewidth',2)
hold on

% Run the simulation and show scattering events. 
Nevents=50;
for ii=1:Nevents
    
    x0=2*(rand-.5)*Max_b;
    y0=2*(rand-.5)*Max_b;
    
    X0=[x0 y0 Z_0 0 0 V_0];
    P.X0=X0;
    P.RK4()
    plot3(P.X(1,:),P.X(2,:),P.X(3,:))
    xlabel('X (m)')
    ylabel('Y (m)')
    zlabel('Z (m)')
end
axis equal



